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Numerical  Simulalion  of  BGK-Bumett  Equations 


ABSTRACT 

Recently’  it  has  been  shown  using  Boltzmann’s  H-Theorem  that  the  conventional  Burnett 
equations  violate  the  second  law  of  thermodynamics,  and  hence  must  not  be  employed  for  fluid 
dynamic  simulations.  To  overcome  this  difficulty,  a  new  set  of  equations,  designated  the  BGK- 
Bumett  equations  was  derived  recently  by  the  authors.  A  second-order  distribution  function  was 
derived  by  employing  the  Chapman-Enskog  expansion  on  the  BGK-Boltzmann  equation. 
Moments  of  the  BGK-Boltzmann  equation  with  the  collision  invariant  vector  using  the  second- 
order  distribution  function  yield  the  BGK-Bumett  equations.  It  has  been  shown  by  the  authors 
that  the  BGK-Bumett  equations  are  stable  to  small  wavelength  disturbances  and  that  they  yield 
results  consistent  with  the  second  law  of  thermodynamics.  In  order  to  prove  that  these  equations 
are  indeed  entropy  consistent,  it  is  shown  that  the  second-order  distribution  function  does  not 
violate  Boltzmann’s  H-Theorem.  This  new  set  of  equations  must  be  used  for  computing 
hypersonic  flows  at  moderate  Knudscn  numbers.  The  BGK-Bumett  equations  are  employed  to 
compute  the  hypersonic  shock  structure.  The  results  of  the  computations  show  that  under  certain 
flow  conditions,  the  conventional  Burnett  equations  violate  the  second  law  of  thermod3mamics 
while  the  BGK-Bumett  equations  provide  entropy  consistent  results. 
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introduction 

In  one  of  the  earliest  attempts  to  solve  the  Burnett  equations,  Fiscko  and  Chapman^  solved  the 
hypersonic  shock  structure  problem  by  relaxing  an  initial  solution  to  steady  state.  They  obtmned 
solutions  for  a  variety  of  Mach  numbers  and  concluded  that  the  Burnett  equations  do  indeed 
describe  the  normal  shock  structure  better  than  the  Navier-Stokes  equations  at  high  Mach 
numbers.  The  equations  were  however  unstable  when  the  grids  were  made  progressively  finer.  In 
a  subsequent  attempt,  Zhong^  showed  that  the  equations  could  be  stabilized  by  adding  a  few  ad 
hoc  super-Burnett  terms  (linear  third  order  terms,  in  order  to  m^tain  second  order  accuracy)  to 
the  stress  and  heat  transfer  terms  in  the  Burnett  equations.  This  set  of  equations  was  termed  the 
“Augmented  Burnett”  equations.  The  Augmented  Burnett  equations  did  not  present  any  stability 
problems  when  they  were  used  to  compute  the  flow  parameters  in  the  hypersonic  shock  structure 
and  hypersonic  blunt  body  problems.  However  attempts  at  computing  the  flowfields  for  blunt 
body  wakes  and  flat  plate  boundary  layers  even  with  the  Augmented  Burnett  equations  have  not 
been  entirely  successful.  It  has  been  conjectured  by  Chapman  et.  al*’*  that  this  instability  may  be 
due  to  the  fact  that  the  Burnett  equations  violate  the  second  law  of  thermodynamics  at  higher 
Knudsen  numbers. 

The  main  objectives  of  the  present  work  are: 

a)  To  formulate  a  methodology  for  deriving  and  integrating  a  new  set  of  entropy 
consistent  Burnett  equations  (designated  as  BGK-Bumett)  that  can  be  extended  to  higher 
dimensions. 

b)  To  check  if  the  constitutive  relations  for  the  BGK-Bumett  stress  and  heat  transfer  terms 
correctly  model  the  flow  properties  at  high  Knudsen  numbers. 
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To  computationally  check  if  the  BGK-Bumett  equations  are  stable  to  small  wavelength 
disturbances. 


The  second-order  distribution  function  is  derived  by  representing  the  collision  integral  in  the 
Bhatnagar-Gross-Krook  (BGK)  form  and  considering  the  first  three  terms  in  the  Chapman- 
Enskog  expansion.  The  BGK  form  of  the  collision  integral  assumes  that  the  collision  processes 
are  predominantly  binary  in  nature.  In  deriving  the  second  order  distribution  function,  an  as  yet 
unanswered  question  is  the  approximation  for  the  material  derivatives  that  appear  in  the  second 
order  terms.  In  a  recent  attempt®,  the  first  and  second  order  distribution  functions  were  obtained 
iteratively  by  perturbation  analysis  of  the  1-D  BGK-Boltzmann  equation.  In  this  analysis  the  Euler 
equations  were  used  to  approximate  the  material  derivatives  in  the  first  order  distribution 
function.  Moments  of  the  BGK-Boltzmann  equation  with  the  collision  invariant  vector  and  the 
first  order  distribution  function  yield  the  Navier-Stokes  equations.  In  order  to  keep  in  step  with 
the  iterative  process  the  Navier-Stokes  equations  were  used  to  approximate  the  material 
derivatives  in  the  second  order  terms.  The  BGK-Bumett  equations  are  obtained  by  taking 
moments  of  the  BGK-Boltzmann  equation  with  the  collision  invariant  vector  and  the  second-order 
distribution  function.  This  set  of  equations  contains  all  the  stress  and  heat  transfer  terms  reported 
by  Fiscko  and  Chapman^  and  has  additional  terms  which  are  similar  to  the  Super-Burnett  terms. 

In  order  to  prove  analytically  that  these  equations  are  indeed  entropy  consistent  it  has  been  shown 
that  the  second-order  distribution  function  does  not  violate  the  H-theorem.  Since  the  definition  of 
the  Maxwellian  and  the  higher-order  distribution  functions  used  in  deriving  the  BGK-Bumett 
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equations  takes  into  account  the  internal  energy  of  the  molecules  and  further  does  not  assume  the 
molecules  to  be  monoatomic.  a  modified  H-function  has  been  formulated  to  prove  the  H-theorem. 
In  order  to  ensure  that  the  entropy  gradient  remains  positive  throughout  the  flow  field,  a  set  of 
boundary  conditions  has  been  derived  using  the  Gibbs  entropy  equation.  It  has  been  shown  that  a 
positive  entropy  production  can  be  ensured  by  setting  the  heat  transfer  terms  to  zero  at  stations 
far  upstream  and  far  downstream  of  the  shock. 


1-D  BGK-BOLTZMANN  EQUATION 

The  1-D  Boltzmann  equation  can  be  written  as  follows,  using  the  BGK  approximation  for  the 
collision  integral  J(fif)- 

^  +  v-^  =  j(f,f)  =  v{f(®^-f)  - 

In  the  above  equation,  f  denotes  the  distribution  function,  v  denotes  the  molecular  velocity,  J(f,f) 
denotes  the  collision  integral,  and  v  denotes  the  collision  frequency.  In  this  representation  the 
non-linear  collision  integral  is  approximated  by  a  single  relaxation  time  model.  This  approximation 
assumes  that  any  non-equilibrium  distribution  function  will  settle  down  to  the  equilibrium 
distribution  exponentially. 


ZEROTH-ORDER  (MAXWELLIAN)  DISTRIBUTION  FUNCTION 

For  the  special  case  of  collision  equilibrium,  the  distribution  function  takes  the  form  shown  in  eq. 

(2).  It  can  been  shown  that  this  is  both  a  necessary  and  sufiicient  condition  for  collision 
equilibrium. 
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Df  denotes  the  number  of  translational  degrees  of  freedom,  u  denotes  the  fluid  velocity  and 

R  = _ —  I  denotes  the  internal  energy  that  accounts  for  the  energy  contribution  due  to  all  the 

2RT‘ 

non-translational  degrees  of  freedom  and  Iq  denotes  the  average  internal  energy,  and  is  given  by 
the  expression 

(2^Dr)-gfRT  (3) 

HIGHER-ORDER  DISTRIBUTION  FUNCTIONS 

The  various  higher-order  distribution  functions  are  obtained  by  representing  them  as  a  Chapman- 
Enskog  asymptotic  series  expansion. 

f=:f(0)+v(‘)+^2f(2)^ . .  (4) 

where,  ^  denotes  the  Knudsen  number.  The  first  two  terms  of  the  Chapman-Enskog 

expansion  give  rise  to  the  first-order  distribution  function  f  =  f^®^  and  the  first  three 

terms  of  the  Chapman-Enskog  expansion  give  rise  to  the  second-order  distribution  function, 
f  =  f  {®)  +  (*)  +  ^^fW .  In  order  to  obtain  exact  analytical  expressions  for  the  first  and  second- 

order  distribution  functions  the  BGK-Boltzmann  equation  is  non-dimensionalized  by  defining  the 
following  non-dimensional  variables 
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X  .  V 

L *  ^  C 
^  '^rms 


^nns 


'rms 


(5) 


where,  L  denotes  the  characteristic  length,  Cnns  denotes  the  root  mean  square  molecular  velocity 
and  Xqo  denotes  the  free  stream  mean  free  path.  On  substituting  the  non-dimensional  variables,  eq. 
(1)  takes  the  following  form: 


ef 

—  +  V- 

dt 


THE  FIRST-ORDER  DISTRIBUTION  FUNCTION 

Substituting  the  first  two  terms  of  the  Chapman-Enskog  expansion,  f  =  f^®^  in  the  BGK- 


Boltzmann  equation  and  equating  like  powers  of  the  Knudsen  number  yields: 


(7) 


On  expressing  f  and  substituting  in  eq.  (7)  the  following  expression  is  obtained.  In 

this  expression  the  1-D  Euler  equations  have  been  used  to  express  the  time  derivatives  in  terms  of 


the  spatial  derivatives. 

1 

II 

A<'>(..)|- 

where. 

a(‘)(I,c)  = 

5c  4Ic(y  - 1 

2p  (3-t) 

Pc^(3-t)  + 

dx 


(8) 


-c 


2  (3-r) 


(9) 


(10) 
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and  c  denotes  the  peculiar  or  thermal  velocity,  c=v-u.  The  first-order  (Kstribution  function  satisfies 
the  property  .  Hence  =  0 . 


SECOND-ORDER  DISTRIBUTION  FUNCTION 

The  second-order  distribution  function  is  obtained  by  considering  the  first  three  terms  in  the 
asymptotic  Chapman-Enskog  expansion  for  the  distribution  function.  Substituting  the  expression 
for  the  second  order  distribution  function  in  the  non-dimensional  1-D  BGK-Boltzmann  equation 
and  equating  like  powers  of  ^  yields: 


_5 

at 


From  the  above  equation  the  following  equation  for  f^^^  is  obtained: 


(11) 


r(^) 


(12) 


Since  the  field  vector  Q  is  the  same,  for  the  Euler  ,  Navier-Stokes  and  BGK-Bumett  equations, 
moments  of  the  distribution  function  with  the  collision  invariant  vector  'F  must  be  the  same  for 
any  distribution  function.  Hence,  the  second-order  distribution  function  must  satisfy  the  property 

('F,  f)  =  ^'F,  f^^^y  This  condition  translates  to  the  following  equation: 


=  0 


(13) 


\  vLat'"  ^  ax 

An  expression  for  the  second-order  distribution  function  satisfyng  the  above  equation  is  given  by: 


f(0)  +5f(O)4,0)  _ir|[f(0)4,(>))+A[uf(0)4>(l))+A(cf(»)3.«) 


(14) 
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where, 

$a).-±x(»(T.c)f-±AW(T.o)|  a^) 


^C  +  — Tc  +  03C^ 

P  3 


(16) 


A^^)(T,c)  =  [3840^  +051+06] 


(17) 


In  the  above  expressions  T  =  —  and  c  denotes  the  peculiar  or  thermal  velocity,  c  -  v-u.  The 

coefficients  0i,  i  =  1,2.. ..6  are  functions  of  the  specific  heat  ratio  Y  •  The  exact  expressions  for 
these  coefiBcients  are  given  in  the  appendbc. 


THE  BGK-BURNETT  EQUATIONS 

The  various  fluid  dynamics  equations  are  obtained  by  taking  moments  of  the  BGK-Boltzmann 


equation  with  the  collision  invariant  vector 


(18) 


00  00 

The  moment  of  f  is  defined  as(T,f)  =  j  J'Ffdvdl.  On  taking  moments  of  the  BGK- 

0  -CO 

Boltzmann  equation,  the  following  generic  equation  is  obtained. 


f) + A(4-.  vf) = ('p.  v(f  0)  -  f)) = 0 


(19) 


On  substituting  the  zeroth,  first  and  second-order  distribution  functions  the  Euler,  Navier-Stokes 
and  BGK-Bumett  equations  are  obtained  respectively. 


^  I  JfTdvdl  J  JvfFdvdl  =v|  -f)dvdl  =  0 


Moments  of  the  BGK-Boltzmann  equation  with  using  the  second-order  distribution  function 
yield  the  BGK-Bumett  Equations. 

vp.fW +5fW  +52f(2)  +yW  - f)  'j  (21) 

The  1-D  BGK-Bumett  equations  are  represented  in  conservation  law  form  as: 

dQ  dG'  dG''  ac®  . 

at  ax  ax  ox 

The  elements  of  the  field  and  flux  vectors  are; 

0 

..^N-S  ^„N-S 

-uXx  +  q 

The  Navier-Stokes  and  BGK-Bumett  stress  and  heat  transfer  terms  are  given  by  the  following 


Q  =  pu  .  G*  =  p  +  pu^  , 
pe  pu  +  pue 
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The  coefficients  and  in  the  expressions  for  the  BGK-Bumett  stress  and  heat 

transfer  terms  are  functions  of  ‘  y  ’  and  are  given  in  the  appendbc.  The  BGK-Bumett  stress  and 

heat  transfer  terms  have  two  sets  of  terms  of  orders  and  .  The  former  results  when  the 
Euler  equations  are  used  to  express  the  material  derivatives  in  terms  of  the  spatial  derivatives  and 
the  latter  is  obtained  when  the  Navier-Stokes  equations  are  used  to  express  the  material 

derivatives  in  terms  of  the  spatial  derivatives.  When  only  terms  of  order  p  are  considered  it  is 
observed  that  the  derivatives  in  the  expressions  for  the  stress  and  heat  transfer  terms  are  identical 
to  the  derivatives  reported  by  Fiscko  and  Chapman^.  The  coefficients  of  these  derivatives  are, 
however,  very  different  to  those  in  Ref  3.  Table  1.  shows  the  comparison  between  the  BGK- 
Bumett  coefficients  and  the  coefficients  of  the  Burnett  equations  in  Ref  3.  When  terms  up  to  the 

order  p^  are  considered  it  is  observed  that  there  are  many  non-linear  terms  (products  of 
derivatives)  in  addition  to  linear  third-order  derivatives.  These  derivatives  are  similar  to  the  super- 
Burnett  derivatives  in  Ref  3.  Table  2.  and  Table.  3.  show  the  comparisons  between  the 
coefficients  of  the  BGK-Bumett  derivatives  and  the  coefficients  of  the  augmented  Burnett  terms 

evaluated  by  Zhong^. 


I 
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LBVEARIZED  STABILITY  ANALYSIS 

It  has  been  shown  by  Bobylev’  that  the  conventional  Burnett  equations  are  not  stable  to  small 
wavelength  disturbances.  Hence  the  conventional  Burnett  equations  tend  to  blow  up  when  the 
mesh  sizes  are  made  progressively  finer.  In  order  to  investigate  the  stability  aspects  of  the  BGK- 
Bumett  equations  a  model  problem  is  considered  which  studies  the  response  of  a  uniform  gas  to  a 
1-D  periodic  perturbation  wave.  The  initial  density,  temperature  and  velocity  of  the  undisturbed 
gas  at  time  t  =  0  are  po,  Tq,  and  Uq  =  0  respectively.  At  t  =  0  the  gas  is  perturbed  such  that: 


( 

ioox^ 

P  =  P0 

1  +  Cie^ 

(27) 

( 

iCDX^ 

H 

II 

1  +  C2e^ 

(28) 

> 

1 

^  iCDX^ 

u  =  -JkTq 

Cje^ 

(29) 

1 

1  / 

Since  the  perturbations  are  assumed  to  be  small  the  magrutudes  of  the  coefficients  in  the 
expressions  (27)-(29)  are  required  to  satisfy  the  inequality  jCk]  « 1,  (k  =  1,2,3). 


The  characteristic  length  Lq  =  — =  0.783X,  where  X  denotes  the  mean  fi-ee  path.  The 

PoV^^T'o 


271  X 

non-dimensional  circular  fi-equency  o  =  - r  =  4.92—  =  4.92Kn .  Introducing  the  perturbed 


quantities  in  the  continuity,  momentum  and  energy  equations  and  simplifying  yields: 


Numerical  Simulatioa  of  BGK-Bumett  Equations 


av  ,,  av-  ,,  a^v.„  a"v,,,  a"v_^ 

at*  ^  ax'  ^  ax'2  ax'^  ax'^ 


■p']  0  1  0 

V'=-  u'  >,  Ml  =  1  0  1  ,  M2 

TJ  0  (y-1)  0_ 


-(3-y)  0  ,  M3  = 


0 

PrJ 


R  R 

0  bW(T  - 1)  0 


0  0 

(9) 

0  ^ 
R 

0  0 


bP)(T-I) 


The  non-dimensional  initial  conditions  for  eq.  (30)  can  be  denoted  in  vector  form  as 

V'l  =  Ve‘“*  ,  where  x'  =  .  Let  us  assume  the  solution  of  the  above  equation  to  be  of  the 

Lq 


V=  Ve'^^’c*'’ 


(31). 


where  t'  =  .  The  complex  variable  4»  =  a  +  iP .  a  denotes  the  attenuation  coefficient  and  P 

Mo 

denotes  the  dispersion  coefficient.  For  stability  a  <  0  as  L  decreases  or  in  other  words  the  flow 
must  attenuate  as  the  Knudsen  number  increases.  Substituting  eq.  (31)  in  eq.  (30)  and  simplifying 
yields  eq.  (32)  when  Euler  equations  are  used  to  express  the  material  derivatives. 

|<1)I  +  icdMi  -  (D  ^M2  -  ico  ^M3  jVoe‘“*  e*^  =  0  (32) 

For  a  non-trivial  solution  the  following  condition  must  be  satisfied 


+  itoMi  -  CO  ^M2  -  ico  ^M3  =  0 
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When  the  Navier-Stokes  equations  are  used  to  express  the  material  derivatives  in  terms  of  the 
spatial  derivatives  the  following  equation  is  obtained: 


|c>I  +  icoMi  —  ©  ^M2  “  i©  ^M3  +  ©  ^M4  jVoe**”*  e*^  —  0  (34) 

For  non-trivial  solutions  the  following  condition  must  be  satisfied: 


0  ^  A 

OI  +  i©Mj  —  ©  M2  “  i©  M3  +  ©  M4 


=  0 


(35) 


The  trajectory  of  the  roots  of  the  characteristic  equations  (33)  and  (35)  is  plotted  on  the  complex 
plane  on  which  the  real  axis  denotes  the  attenuation  coeflScient  and  the  imaginary  axis  denotes  the 
dispersion  coefficient.  For  stability  it  is  required  that  the  roots  lie  to  the  left  of  the  imaginary  axis 
as  the  Knudsen  number  increases.  Fig(s)  1-4  show  the  trajectory  of  the  roots  of  the  characteristic 
equation  as  the  Knudsen  number  increases.  From  the  plots  it  is  observed  that  unconditional 
stability  is  guaranteed  only  when  the  Navier-Stokes  equations  are  used  to  express  the  material 
derivatives  in  terms  of  the  spatial  derivatives.  It  must  be  noted,  however,  that  the  linear  stability 
analysis  does  not  consider  the  many  non-linear  terms  -  powers  and  products  of  derivatives  -  that 
are  present  in  the  BGK-Bumett  stress  and  heat  transfer  terms.  Hence,  this  analysis,  is  at  best  only 
a  necessary  condition  for  the  stability  of  these  equations.  A  more  rigorous  proof  of  the  stability  of 
these  equations  involves  verifying  the  Boltzmann’s  H-Theorem. 


BOLTZMANN’S  H-THEOREM 

The  BGK-Bumett  equations  must  satisfy  the  second-law  of  thermodynamics.  There,  however,  is 
no  acceptable  definition  of  entropy  for  a  gas  in  a  state  of  non-equilibrium.  Physical  intuition  tells 
us  that  an  isolated  system  will  evolve  fi'om  an  arbitrary  initial  state  to  a  state  of  equilibrium. 
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Numerical  Simulation  of  BGK-Biimett  Equations 


Boltzmann’s  H  theorem  formalizes  this  notion,  and  also  makes  explicit  the  manner  in  which  this 
evolution  proceeds.  A  spatially  homogenous  gas  is  defined  as  one  in  which  the  density  does  not 
vary  with  position.  Boltzmann’s  H  theorem  states  that  for  a  spatially  homogenous  gas  the 

inequality,  — (H)^0,  must  be  satisfied  when  the  gas  approaches  equilibrium.  The  quantity  H 
dt 

00 

which  is  shown  to  be  the  kinetic  theory  equivalent  of  entropy’  is  defined  asH=jflnfdv.In 

-CO 

arriving  at  this  definition  of  the  H-function  Boltzmann  made  the  following  assumptions: 

a)  The  molecules  comprising  the  gas  do  not  have  any  internal  energy.  Hence  the  H 
function  was  defined  only  over  the  range  of  molecular  velocities. 

b)  The  gas  was  assumed  to  be  monoatomic. 

Since  our  definition  of  the  Maxwellian  and  the  first  and  second  order  distribution  functions  takes 
into  account  the  energy  contribution  due  to  the  various  non-translational  degrees  of  fi'eedom  and 
further  does  not  assume  the  gas  to  be  monoatomic,  the  definition  of  H  must  be  modified  to 
account  for  these  differences.  The  modified  definition  of  the  H  function  can  be  shown  to  reduce  to 
the  classical  (Boltzmann)  definition  of  H  for  the  specific  case  of  a  monoatomic  gas. 

MODIFIED  H-FUNCTION 

The  change  in  entropy  in  classical  thermodynamics  is  given  by  the  following  expression: 
(s2-Si)=Cvln^-Rln^ 

^  Ti  Pi 


'  It  must  be  noted  that  entropy  according  to  classical  thermodynamics  is  defined  only  for  equilibrium  systems.  The 
quantity  ‘H’ ,  however,  is  defined  even  for  non-equilibrium  systems. 
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Numerical  Simulation  of  BGK»Bumett  Equations 


where  the  subscripts  *1’  and  ‘2’  denote  thermodynamic  variables  at  equilibrium  stations  far 
upstream  and  far  downstream  respectively.  The  absolute  entropy  is  ^ven  by  the  following 
expression  up  to  an  additive  constant. 

s=CvlnT-Rlnp  +  R^  (3'7) 

The  above  expre.ssion  can  be  cast  in  the  form 


s=-R 


I  V9 

lnp  +  7 - 

(y  - 1) 


(38) 


where  B  =  — ^  We  now  need  to  devise  a  method  to  arrive  at  the  above  expression  from  the 
2RT‘ 

Maxwellian  distribution  function  which  is  also  the  equilibrium  distribution  function.  The  1-D 
Maxwellian  distribution  function  is  given  by  the  expression 


f(o)  =  F = -T 

Iq  V  rt  |_  iQ 


(39) 


^3  —  y) 

where  the  average  internal  energy  In  =  ^  •  Equation  (40)  can  be  rewritten  as: 

4B(y  -  1) 


lnf^®^  =  lnF  = 


+  — Inp  4  In^ 


<Y-0 


+  (2|3u)v 


pu-* 


(41) 


On  rearranging  the  terms  in  eq.  (41) 


lnf(°)  -231^^-^  = 
(3-y) 


Inp  +  -InP  +  Inf  ^  ^L|  -  pu^ 

2  l(3-y)VnJ 


-2P 


^  2^ 
V  2/ 


(42) 


(2pu)v 
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Equation  (42)  is  indeed  a  linear  combination^  of  the  collision  invariants.  Hence, 


(3-y)JI  a  &J/ 


(3-y)  \  dt  dx  , 


The  above  equation  can  be  recast  in  the  following  form  by  making  use  of  the  identity 


af(o)  5f^®)"l\  ^ 

- +  v -  )=  0. 

dt  dx  I 


(3-y)A  ^  ^  J/ 

yl  P-T)jl  a  aj 

0  — OO 


On  simplifying  ,  eq.  (44)  can  be  expressed  as 


J  J(|.vD|’f<'’)lnf<»)*Mf<'’)lnpldvdI  -  0 


The  above  equation  can  be  expressed  in  the  following  compact  form 


i(„(o))4(h(o))=„ 


where  the  functionals  and  are  defined  in  (47)  and  (48)  as: 


h(°)=J  j|^f<°hnf^°)+^^^f^°hnp  dvdl 


^  It  can  be  shown  that  the  moments  of  the  BGK-Boltzmann  equation  with  any  linear  combination  of  the  collision 
invariants  equals  zero.  Hence  the  need  to  express  equation  (36)  as  a  linear  combination  of  collision  invariants. 
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HV 


(0) 


=  Jj|f(»)|„f(«).^f<»)|nP 

A  ^ 


0  -CO 


dv  dl 


(48) 


On  evaluating  the  moments  in  eq.  (47)  the  following  expression  is  obtained  for  H 


(0) 


H' 


(«)=, 


1  .9 


(49) 


On  comparing  eq(s).  (49)  and  (38), 


ps 


=  -rh(°^ 


(50) 


The  additive  constant  Co  is  a  function  of  the  specific  heat  ratio  y  and  the  gas  constant  R.  The 
definition  of  the  H-fiinction  can  be  extended  to  any  distribution  function.  Accordingly  the 
functionals  H  and  Hy  are  defined  as  shown  in  eq(s).  (51)  and  (52). 


H  =  |J 


0  -00 


flnf  +  — ^flnP 


dvdl 


(51) 


00  00  p 


flnf  +  l^-^flnp 
2(r-i) 


dv  dl 


(52) 


00  00 


It  can  be  seen  that  the  expression  for  H  reduces  to  the  classical  definition  H  -  a  flnf  dv  dl,  for 


0  —00 


the  specific  case  of  a  monoatomic  gas  (i.t.  y  »  )^  )•  For  a  spatially  inhomogenous  gas  Grad*  has 
shown  that  the  following  inequality  must  be  satisfied  when  the  gas  approaches  equilibrium. 


*  This  relation  establishes  a  link  between  Boltzmann’s  H-Theorem  and  the  cl^ical  thermodjTiamics  concept  of 
entropy.  It  must  be  pointed  out  that  there  is  no  rigorous  justification  to  extend  this  definition  to  include  H-fimctions 
derived  fi-om  higher-order  distribution  functions. 
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(53) 


The  Hbove  inequality  will  be  shown  to  be  true  for  the  first  snd  second-order  distribution  functions. 
The  first  and  second  order  distribution  functions  are  given  by  expressions  (54)  and  (55). 


f  =  f(°)  +  ^f(°)(D(^)  =  f(°)(l  + 

(54) 

f  =  f(0)  +  ^f(®)(I)(0  +  ^2f(0)(j,(2)  _  |.(o)^j  ^  ^,j,(l)  ^  ^2(j,(2)J 

(55) 

On  expressing  ln(f)  as  a  Taylor  series  and  considering  only  terms  up  to  the  second  power  in  the 
Knudsen  number  the  following  approximations  are  obtained  for  the  first  and  second  order 
distribution  functions. 


•'))  j  =  InfW  +  - 

{')  +  4’<1>W)1  =  In  H")  +  (i<l>W  + 1 


(56) 


(57) 


The  H-Balancc  Equation  for  the  Boltzmann  equation  is  obtained  by  evaluating  the  moments  in  the 
following  equation. 


00  00 


(3-y)  . 


dv  dl 


(58) 


0  -OO 

On  substituting  eq(s).  (54)  and  (56)  in  the  above  equation  and  retaining  only  terms  up  to  the  first 
power  in  ^  (Knudsen  number)  the  right  hand  side  (RHS)  of  the  above  equation  takes  the  form; 


^Jiimeiicai  Simuladon  of  BGK-Bmuett  Equations 


00  00 


j  j 


0  —00 


(3-r)  . 


dv  dl 


(59) 


The  above  integral  equals  zero  as  moments  of  the  first-order  distribution  function  with  a  linear 


combination^  of  the  collision  invariants  are  being  evaluated.  Since  the  first-order  distribution 
function  satisfies  the  property  =  eq.  (54)  equals  zero.  Hence  the  first  -order 

distribution  function  satisfies  Boltzmann’s  H-Theorem.  On  evaluating  eq.  (58)  up  to  the  second 
power  in  ^ ,  the  RHS  of  the  equation  takes  the  following  form: 


J  J  lnf(®)  -  ^f(®)<l>W ^  J  I  dv  dl 

0-<o  ^ 


(60) 


0  •CO 


The  second  integral  equals  zero  as  the  second-order  distribution  function  satisfies  the  moment 
property  ^'F;  =  ^'F;  =  0 .  For  reasons  given  earlier  moments  of  the  second- 


order  term  with  any  linear  combination  of  the  collision  invariants  equals  zero.  Hence  eq.  (60) 
simplifies  to  the  following  : 


-J  dvdl 


(61) 


0  —00 


On  substituting  for  fi^om  eq.  (8)  and  evaluating  the  moments 


*  Moments  of  the  first  and  second-order  distribution  function  with  any  linear  combination  of  the  collision 
invariants  equals  zero. 
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It  must  be  noted  that  eq.  (57)  is  always  less  than  zero.  Hence  the  second-order  distribution 
function  satisfies  Boltzmann’s  H-Theorem.  It  has  been  proved  conclusively  that  both  the  first  and 

i 

second-order  distribution  functions  satisfy  the  inequality 

|(H).|;(Hv).0 

As  a  consequence  the  BGK-Bumett  equations  which  are  obtained  fi-om  the  BGK-Boltzmann 
equation  by  using  the  second-order  distribution  function  are  entropy  consistent ! 


NUMERICAL  EXPERIMENTS 

The  BGK-Bumett  equations  were  numerically  integrated  using  a  hybrid  algorithm  for  the 
hypersonic  shock  structure  problem.  In  the  hybrid  algorithm  the  inviscid  fluxes  were  split  using 
the  KWPS  scheme  and  the  viscous  and  BGK-Bumett  fluxes  were  central  differenced.  The 
objective  of  this  experiment  was  to  test  the  computational  stability  of  the  entropy  consistent 
BGK-Bumett  equations  by  integrating  them  numerically  on  progressively  finer  grids.  In  order  to 
test  the  stability  of  the  algorithm,  the  scheme  was  applied  initially  to  a  coarse  mesh  of  101  grid 
points.  The  number  of  grid  points  was  increased  to  501.  The  reference  parameters  used  for  Argon 
are  similar  to  those  used  by  Zhong*.  The  results  of  the  computations  are  shown  in  Fig(s)  5-12. 
Fig(s)  5-8  compare  the  entropy  plots  of  the  BGK-Bumett  equations  with  the  entropy  plots  of  the 
Burnett  equations  of  Fiseko  and  Chapman’.  Fig(s)  9-12  compare  the  normalized  temperature  and 
specific  entropy  profiles  of  the  BGK-Bumett  and  Navier-Stokes  equations. 
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CONCLUSIONS 

An  entropy  consistent  set  of  BGK-Bumett  equations  has  been  derived  from  first  principles.  These 
equations  have  been  numerically  integrated  to  compute  the  hypersonic  shock  stmcture.  The 
equations  are  computationally  stable  for  the  range  of  grid  points  and  Knudsen  numbers  for  which 
results  are  presented. 


PROPOSED  RESEARCH 

An  attempt  is  being  made  to  compare  the  expressions  for  entropy  balance  obtained  from  kinetic 
theory  (Boltzmann’s  H-Theorem)  and  classical  thermodynamics,  and  account  for  the  differences 
between  the  two  expressions.  Since  the  BGK-Bumett  equations  in  1-D  have  been  proven  to  be 
entropy  consistent  a  similar  analysis  for  2-D  BGK-Bumett  equations  will  be  earned  out. 
Numerical  solutions  for  hypersonic  flow  past  blunt  bodies  will  be  computed  using  the  2-D  BGK- 
Bumett  equations. 
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Burnett  and  BGK-Bumett  Coefficients 


Huler  equations  have  been  used  in  the  BGK-Burnett  equations  to  express  the  tune  derivatives  in 
terms  of  the  spatial  derivatives. 

Table;  1 


Stress  and  Heat 

Transfer 

CoefiBcients 

Fiscko  & 

Chapman  (Air) 

Y  =  1.4 

BGK-Bumett 

(Air) 

Y  =  1.4 

Fiscko  & 

Chapman 

(Argon) 

Y  =  1.6666 

BGK-Bumett 

(Argon) 

’ 

Y  =  1-6666 

1.749 

0.96 

1.749 

0.446 

.P) 

-388.024 

-459.2 

-281.216 

-277.472 

.P) 

388.024 

459.2 

281.216 

277.472 

-257.726 

-5625.0 

-186.784 

-1216 

403.522 

-51660 

292.448 

-938.11 

a(6) 

74.62 

i 

1 

-5625 

54.08 

-1216.0 

10.831 

-21.633 

10.831 

-9.896 

bP) 

-2.269 

0.183 

-2.269 

-0.194 

b(’) 

-2.06 

0.533 

-2.06 

-0.443 

23 


RBK-RKA 


Numerical  of  BGK^Burnett  Equations 


Augmented  Burnett  and  BGK-Bumett  Stress  Coefficients 

Navier-Stokes  equations  have  been  used  to  express  the  time  derivatives  in  terms  of  the  spatial 
derivatives. 

Table:  2 


Stress 

CoefiBcients 

Zhong  (Air) 

Y  =  1.4 

BGK-Bumett 
(Air) 
y  =  1.4 

Zhong  (Argon) 
y  =  1.6666 

BGK-Bumett 
(Argon) 
y  =  1.6666 

a« 

1.749 

0.96 

1.749 

0.446 

aW 

-388.024 

-459.2 

-281.216 

-277.472 

a(’) 

388.024 

459.2 

281.216 

277.472 

a(‘) 

-257.726 

-5625 

-186.784 

-1216 

a(’) 

403.522 

-5166 

292.448 

-938.11 

.w 

74.62 

-5625 

54.08 

-1216 

0 

-459.2 

0 

-277.472 

- - - 

a(*) 

0 

642.88 

0 

462.268. 

a(’) 

63.778 

459.2 

46.222 

277.472 

a(">) 

0 

0.64 

0 

0.888 
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Augmented  Burnett  and  BGK-Bumett  Heat  Transfer  Coefficients 

Navier-Stokes  equations  have  been  used  to  express  the  time  derivatives  in  terms  of  the  spatial 
derivatives. 


Table;  3 


Heat  TransfW 
CoefiBcients 

Zhong  (Air) 

BGK-Bumett 

(Air) 

Zhong  (At) 

BGK-Bumett 

(At) _ 

b(>) 

10.831 

-21.633 

10.831 

-9.896 

b(^) 

-2.269 

0.183 

-2.269 

-0.194 

-2.06 

0.533 

-2.06 

-0.443 

bW 

0 

4.4 

0 

4.666 

b(’) 

0 

1406 

0 

866.84 

b(^) 

0 

-1406 

0 

-866.84 

-2186 

1406 

-2041 

866.84 

b(«)  1 

0 

0.7 

0 

0.833 

b(’) 

0 

-14 

0 

-1.666 

b(>«) 

0 

-22500 

0 

-13870 

bO') 

-179.375 

0 

-130 

0 
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express  the  material  derivatives  in  terms  of  the  spatial  derivatives. 


Numerical  Simulation  of  BGK-Bumett  Equations 


equations  have  been  used  to  express  the  material  derivatives  in  terms  of  the  spatial  derivatives 


Numerical  Simulation  of  P.GK-Bumett  Equations 


Attenuation  Coefficient 


Numerical  Simulation  of  BGK-Bumett  Equations 


used  to  express  the  material  derivatives  in  terms  of  the  spatial  derivatives 


Numerical  Simulation  of  BGK-Bumett  Equations 
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Euler  equations  have  been  used  to  express  the  material  derivatives 
spatial  derivatives. 
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Numerical  Simulation  of  BGK-BumettEquatiiOiis 


Fig.  7  Entropy  variation  across  a  Mach  35  normal  shock  (Argon).  Kn^  =  0.0005 

Euler  equations  have  been  used  to  express  the  material  derivatives  in  terms  of  the 
spatial  derivatives. 
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Fig.  9  Temperature  variation  across  a  Mach  35  normal  shock  (Argon).  Kn^o  =  0.01 

Navier-Stokes  equations  have  been  used  to  express  the  material  derivatives  in  terms  of 
the  spatial  derivatives. 
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Fig.  10  Specific  entropy  variation  across  a  Mach  35  normal  shock  CArgonj 

Navier-Stokes  equations  have  been  used  to  express  the  material  derivatives 
the  spatial  derivatives 


Numerical  Simulation  of  BGK-Bumett  Equations 
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Fig.  11  Temperature  variation  across  a  Mach  35  normal  shock  (Argon).  J^oo  —  u.uz 

Navier-Stokes  equations  have  been  used  to  express  the  material  derivatives  in  terms  of 
the  spatial  derivatives. 
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Appendix:  One 

The  expressions  for  the  coefficients  in  the  BGK-Bumett  flux  vector  are  given  below; 


01  =  23  +  toi,  02=0)2-1.  93  =  -(1  +  “i)*  e4  =  (3-Y)  +  “3 


05  =  -[(y- 0  +  0)3].  06  =  [^^^)+“3 

n,  =  .  02  =  (D(e,  +  62  +  fe,) .  03 = 


04  =  w  ^  ~  +  (y  +  3)05  +  2706) 

8(y-0 


The  coefficients  of  the  stress  and  heat  transfer  terms  are  given  below; 


=  2(2  -  7)01,  =  -2niR,  a^O  =  2niR.  a^'*^  =  -4(^y 


a(^)=-+l! 


+  n2jR,  =  -2O1R,  a^O  =  2R7O1,  a^O  =  2O1R,  a^^®^  -  2(7  l)Oi 
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gimiilation  of  BGK-Bumett  Equations 


=  804  -  403(27  - 1)  -  4O2  -  2O1,  =  O4  -  03(7  - 1) ,  b(^)  =  4O4  -  2O1 

bW  =  2(Oi  +  403(7  -  1)) ,  =  47O3R,  b(^)  =  -47O3R,  b(^)  =  47O3R 

b^*^  =  203(7  - 1) ,  b^^^  =  -403(7  - 1) ,  =  -647O3R 
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Abstract 


In  this  paper  two  different  forms  of  Burnett  equations  are  studied  which  have  been 
designated  as  ‘Augmented  Burnett  Equations’  and  ‘BGK-Bumett  Equations’.  The 
augmented  Burnett  equations  were  developed  by  Zhong  to  stabilize  the  solution  of  the 
conventional  Burnett  equations  which  were  derived  in  1935  by  Burnett  from  the 
Boltzmann  equation  using  the  second-order  Chapman-Enskog  expansion.  In  this 
formulation,  ITie  conventional  Burnett  equations  are  augmented  by  addmg  cd  hoc  third- 
order  derivatives  to  stress  and  heat  transfer  terms  so  that  the  augmented  equations  are 
stable  to  small  wavelength  disturbances.  The  BGK-Bumett  equations  have  been  recently 
derived  by  Agarwal  and  Balakrishnan  from  the  Boltzmann  equation  using  the  Bhatnagar- 
Gross-Krook  (BGK)  approximation  for  the  collision  integral.  These  equations  have  been 
shown  to  be  entropy  consistent  and  satisfy  the  Boltzmann  H-Theorem  in  contrast  to  the 
conventional  Burnett  equations  which  violate  the  second  law  of  thermodynamics.  In  this 
paper,  both  sets  of  Burnett  equations  are  applied  to  compute  a  2-D  hypersonic  flow  over  a 
circular  cylinder  at  Knudsen  numbers  0.001  to  0.1.  The  radius  of  the  cylinder,  which  is  the 
characteristic  length  of  the  body,  determines  the  Knudsen  number.  The  Steger-Warming 
flux-vector  splitting  scheme  is  applied  to  the  convective  inviscid  flux  terms.  Stress  and 
heat  transfer  terms  are  simply  second-order  central-differenced.  Comparison  is  made 
between  the  augmented  and  BGK-Bumett  equations  solutions  and  with  the  Navier-Stokes 
calculations.  Comparison  of  the  solutions  from  the  augmented  Burnett  equations  with  the 
Navier-Stokes  solutions  shows  that  the  difference  is  significant  at  high  Knudsen  number 
(A>i=0.1).  The  solutions  from  the  BGK-Bumett  equations  are  matched  well  with  those 
from  the  Navier-Stokes  equations  and  the  augmented  Burnett  equations  at  lower  Knudsen 
numbers  (K>f=0.001,  0.01).  BGK-Bumett  solutions  are  currendy  underway  at  higher 
Knudsen  numbers. 
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Nomenclature 


et 
Kn 
M 
Pr 
P 
<li 
Re 
R 
T 

7; 

t 

u,  V 

Greek  symbols 
a*  3* 

y^ 

5, 

e, 

K 

n 

p 

y 


Superscript 

third  order  terms  in  augmented  Burnett  equations 
third  order  terms  in  BGK  Burnett  equations 


Subscripts 

jc,  y  derivatives  in  x  and  y  directions 

„  free  stream  quantities 


coefiGcients  of  stress  terms  in  Burnett  equations 

coefiBcients  of  heat  flux  terms  in  Burnett  equations 

coefiBcients  of  stress  terms  in  Navier-Stokes  equations 

coefiBcients  of  third  order  terms  in  BGK  Burnett  equations 

thermal  conductivity 

coeflBcient  of  viscosity 

density 

stress  tensor 

specific  heat  ratio 


total  energy 
Knudsen  number 
Mach  number 
Prandtl  number 
pressure 
heat  flux 

Reynolds  number 
gas  constant 
temperature 
wall  temperature 
time 

velocity  components  in  x,  y,  z  direction 


1.  Introduction 


In  one  of  the  first  attempts  to  solve  the  Bumett  equations,  Fiscko  and  Chapman^^^  solved 
the  hypersonic  shock  structure  problem  for  a  variety  of  Mach  numbers  and  concluded  that 
the  Bumett  equations  describe  the  normal  shock  stmcture  better  than  the  Navier-Stokes 
equations  at  high  Mach  numbers.  However,  in  their  numerical  solution,  they  experienced 
stability  problems  on  finer  grids.  The  linearized  Bumett  equations  were  found  to  be 
unstable  to  small  wavelength  disturbances.  In  a  subsequent  attempt,  Zhong®  stabilized  the 
Bumett  equations  by  adding  a  few  linear  third  order  terms  on  an  otf  hoc  basis.  This  set  of 
equations  was  termed  the  augmented  Bumett  equations.  The  augmented  Bumett 
equations  did  not  present  stability  problems  when  they  were  tqiplied  to  the  hypersonic 
shock  stmcture  and  hypersonic  blunt  body  problems.  However  the  augmented  Bumett 
equations  were  not  entirely  successful  to  compute  the  flowfields  for  blunt  body  wakes  and 
flat  plate  boundary  layers.  Comeaux  et.  al®  have  noted  that  the  linear  stability  analysis  is 
not  sufficient  to  explain  the  instability  of  the  Bumett  equations  with  increasing  Knudsen 
numbers  because  of  many  non-linear  terms  present  in  the  Bumett  equations.  They  have 
conjectured  that  this  instability  may  be  due  to  the  fact  that  the  Bumett  equations  violate 
the  second  law  of  thermodynamics  at  high  Knudsen  numbers. 

The  highly  non-linear  nature  of  the  collision  integral  in  the  Boltzmann  equation  is 
simplified  by  representing  the  collision  integral  in  the  Bhatnagar-Gross-Krook  (BGK) 
form.  Balakrishnan  and  Agarwal®  have  formulated  a  new  set  of  entropy  consistent  1-D 
Bumett  equations  fi'om  the  BGK-Boltzmann  equation  and  using  the  Na^der-Stokes 
equations  to  approximate  the  material  derivatives  in  the  second  order  terms  in  the 
Chapman-Enskog  expansion.  The  material  derivatives  are  thus  expressed  in  terms  of 
spatial  derivatives  using  the  Navier-Stokes  equations.  This  set  of  BGK-Bumett  equations 
contains  all  the  stress  and  heat  flux  terms  reported  by  Fiscko  and  Chapman^’^  and  has 
additional  terms  which  are  similar  to  the  super  Bumett  terms.  Recently,  Balakrishnan  and 
Agarwal®  have  extended  the  1-D  BGK-Bumett  equations  to  2-D  BGK-Bumett 
equations.  In  this  paper,  the  augmented  Bumett  equations®  and  the  2-D  BGK-Bumett 
equations^*^  have  been  used  to  compute  and  compare  the  shock  stmcture  and  other  flow 
properties  for  hypersonic  flow  over  a  blunt  body  in  continuum-transition  regime. 
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2.  Governing  Equations 

The  governing  equations  for  2-D  compressible  viscous  flow  in  Cartesian  coor(flnates  are 


^  ^  ^ 
dt  ^  8x^  dy 


where 


(2) 


E  and  F  are  the  flux  vectors  of  the  conserved  variables  Q  in  the  x  and  y  directions. 
These  flux  vectors  can  be  written  as 


E  =  Ej  +Ey, 

F  =  Fy,  ^3^ 

where  Ej  and  Fj  are  the  invdscid  flux  terms  and  Ey  and  Fy  are  the  viscous  flux  terms  given 
as  follows: 


pu 

0 

pu^+p 

On 

.  ^  = 

ptdv 

Ol2 

(e,+p)u_ 

o„i/+a,2V+^,_ 

pv 

0 

puv 

02, 

pv^+p 

»  ~ 

O22 

(e,  +p)v 

a2,i/+a22V+^j_ 

(4) 


(5) 


In  Eqs.  (4)  and  (5),  the  stress  tensors  and  heat  flux  terms,  o^.  and  q,  are  expressed  as 
follows: 


c,  =  a<»>  +  o;>  +  o"’  +o«’  +•■  ■ + o'”  +  0(A>.“' ), 
?,=?;'”+  9r’+9!”+  9!”  +  - + ^”  +  ). 
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The  zeroth  order  approximation  (n=0)  results  in  the  Euler  equations. 


=  0. 

and  (V 

gr  =  0. 

The  first  order  approjdmation  represents  the  Navier-Stokes  equations.  The  stress  tensors 
and  the  heat  flux  terms  (/i=l)  are  given  as, 

a{V  =  -p(6,«,+52V^), 

=-li(V,+S2i/x). 

and  (^) 

where  (  ),  =  didx  and  (  )^  =  didy.  The  coefficients,  5,  and  Sj  are  given  in  Table  1  for  the 
augmented  Burnett  equations^^)  and  the  BGK-Bumett  equations^. 


■ 

Aug.  Burnett  Eqns. 

BGK-Bumett  Eqns. 
Y=  1.4 

BGK-Bumett  Eqns. 

Y  =1.666 

El 

1.333 

1.6 

1.333 

-0.666 

-0.4 

-0.666 

Table  1.  The  coefficients  in  the  Navier-Stokes  Eqns.  stress  tensors. 

Similarly,  the  second  order  approximation  represents  the  Burnett  equations.  The 
expression  for  stress  and  heat  flux  terms  («=2)  are. 


=  -^(  a^ul  +ajW,v  +051/^ 

P 

RT  RT 

+a,/?r^  H-a,— p„  +a,o— 

RT  ^  R^  ^  RT ^2 

2  Px  "*’®^12  TgPx  '*'^13  'P^*  _2  P^” 

P  P  i  P 

R  rp  R  •p2  \ 

■*^15  TyPy  ^leipTy  ), 


(9) 


5 


=-^(y,^r,i/,+Yj-^r,v^+Y3tt„+Y4«^+y5V^ 
p  r  7 

+76  +77^^“/  +7|“Px“x  +79'“Px^>- 

1  1  N 

+7,o-PA+Yn-PA  X 

P  P  (12) 

and 

=  iL(  y  J 1  r  V,  +  y  ^  1  +  y  +  74V„  +  7  J«:^ 

p  r  7 

+76y^"/  +77y^'’x  +7i”P>.'',  +79“Pj-"x 

1  1  ^ 

+7,o-PxVx+7n-p,",)- 

P  P  (13) 

Both  augmented  Burnett  and  BGK-Bumett  equations  have  same  forms  of  the  stress 
tensor  and  heat  flux  terms.  However  the  two  sets  of  equations  have  different  values  of  the 
coeflBcients.  The  coefficients  are  compared  in  Table  2. 
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Aug.  Burnett  Eqns. 
Hard-Sphere  gas 

BGK-Bumett  Eqns. 
Y=1-4 

BGK-Bumett  Eqns. 

Y  =1.666 

a, 

1.199 

-2.24 

-2.222 

ot. 

0.153 

-0.48 

-4.444 

ot, 

-0.600 

0.56 

1.111 

-0.115 

-1.20 

-0.667 

a. 

1.295 

0.0 

0.0 

-0.733 

0.0 

0.0  . 

a. 

0.260 

-19.6 

-5.833 

a. 

-0.130 

-5.6 

0.0 

Cto 

-1.352 

-1.6 

-1.333 

a,n 

0.676 

0.4 

0.667 

ot,, 

1.352 

1.6 

1.333 

a„ 

-0.898 

-19.6 

-5.833 

Ct,i 

0.600 

-18.0 

-4.5 

a,^ 

-0.676 

-0.4 

-0.667 

a„ 

0.449 

-5.6 

0.0 

a„ 

-0.300 

-6.0 

-0.667 

Pi 

-0.115 

-1.4 

-1.667 

P7 

1.913 

-1.4 

-1.667 

Pi 

0.390 

0.0 

.  CLO 

P- 

-2.028 

-2.0 

-2.0 

Pi 

-0.900 

2.0 

2.0 

P« 

2.028 

2.0 

2.0 

P7 

-0.676 

0.0 

0.0 

Y, 

10.830 

-25.241 

-11.101 

Y, 

0.407 

-0.2 

-1.0 

Y, 

-2.269 

-1.071 

-1.384 

Y^ 

1.209 

-2.0 

-2.0 

Y, 

-3.478 

-2.8 

-3.333 

Yk 

-0.611 

-7.5 

-6.5 

Y, 

11.033 

-11.0 

-5.667 

Y. 

-2.060 

-1.271 

-1.051 

Yo 

1.030 

1.0 

1.0 

Y.r, 

-1.545 

-3.0 

-3.0 

Yu  .. 

-1.545 

-3.0 

-3.0 

Table  2.  The  coefficients  in  the  Burnett  and  BGK-Bumett 
Eqns.  stress  tensor  and  heat  flux. 


The  third  order  approximation  (w=3)  represents  the  super  Burnett  equations.  However, 
not  all  of  the  third  order  terms  of  the  super  Burnett  equations  are  used  in  the  augmented 
Burnett  and  the  BGK-Bumett  equations.  In  the  augmented  Burnett  equations,  the  third 
order  terms  are  employed  ad  hoc  basis  to  obtain  stable  numerical  solutions  while 
fnaintaining  second  Order  accuracy  of  the  solutions.  The  third  order  terms  in  the 
augmented  Burnett  equations^^)  are  given  as, 

„3 

Cfjj^  ~  )> 

P  (14) 

m3 

an  V  +a„tt^  +«,,«„  ), 

P  (15) 


~  P(.  y  llTxxx  "^y  Xt^xyy  Y 13  Y 13  "“Pw  )» 

PP  P  P  (17) 

and 

qM  =ii-/?(y,,r^+y,,r^+y,3-^p^+y,3-^p^  ). 

PP  P  P  (18) 

The  superscript  '{a)'  denotes  augmented  Burnett  terms.  The  coefficients  in  stress  and 
heat  flux  terms  are  given  in  Table  3. 


0.11 

0.2222 

a„ 

-0.1111 

3. 

0.1667 

Yu 

0.6875 

Y„ 

-0.625 

Table  3.  The  coefficients  in  the  augmented  Burnett  Eqns. 

The  BGK-Bumett  equations  have  more  additional  third  order  terms  than  the 
augmented  Burnett  equations.  These  are  not  added  on  an  hoc  basis  but  are  derived 
from  the  third  order  Chapman-Enskog  expansion  of  the  BGK-Boltzmann  equation.  The 
third  order  terms  in  the  BGK-Bumett  equations^*)  are  given  as 


RT  . 
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+07pyV„ +e,p^ii^ +04P^v^ ) 

„3 

+-^(e9“x  +3e,oi/"vy +e„M.vj  -Q^u^ui  -2Q^u,u,v, 

p 

-04i/,v^  +e,oV^  -e,2V^z/J  -20,2«A''y  ) 

H  ^XlPx^xx  ‘^^n^x'^yy  '^^\4^yTxx  "^^U^y^yy  )> 

P 

og>  =-^«r(e.v^+9,v^+e,B„+e,i<„ ) 

L|3 

-^— (O.PyV^  +05P,.«^  +66Py>'« 

+  07P^I/^  +0,P,%  +04P*«»  ) 

^3 

+^(09^^  +30,ovjtt,  +0„v^i/*  -04V^v2  -204V^v,j/^ 
-04V^2/^  +0,oZ/^  -012“xV^  -20,2V,tt,t/:.  -Qu^X  ) 
+^/?(0,3V^r^  +0,3Vy7’„  +0i4«x^3y  +0i4«x^x*  ). 

0I5“x9-  +«>»-  +®IJ  V  +  ''»x  ) 

P 

- (  06Py“x*  +0.6P,V,,  +PyU^ 

P  P 

+  PxVx.  +0l6Px“v  +®6PxV^  ) 

+v,  Xe4"x  +26,2WxV^  +207l/^V, 

+  0X+07V*  +04^^  ) 

+  '^5-i?(  0j7U^7^  +0,71/^?^  ■*'0|7'^x^  '^017'^x^  )» 

+e,.r^-e„ip,r„  -e,.ip.7; ) 

PP  P  P 

n  ^ 

+ - (  0,9W,U„  +02o«xV^  +06«x*^3y  +021  V« 

/>p 


(19) 


(20) 


(21) 
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and 


“*"®23^y^xc 

“^®6^y^>y  “^®23^x^xc  “^®24^*^ay  ) 

-■^(-P,  +\t.  X^y.  +2eMV,  +2e„i<,v. 

PP  P  T 

+  0j7tfy  +0I7Vjj  +0i3Vj,  ) 

+— |(9..7;r,+0„r.r„). 

ppT 

,<"  =iilK(e..r  +0„r^ -K)-pyT„ -0,  ip,?^ ) 

PP  P  P 

n  ^ 

+— (  ei9V,,V  +  02oV^«^  +  06V^V„  +  02,tf,V^ 

PP 

+  e^U^U^  +07tf,V„  +023V,«^  +624V,V^ 

+  06V,U„  +023l/7,«^  +e24«yV^  +06","»  ) 

-■'^(-P,  +^3;  X0„x’  +20„..v,+20„u,v. 

/^P  P  ^ 

+  0,7!/^,  +017^,  +013^^  ) 


+“( 9, .j;n. +«..?',';)• 


(22) 


(23) 


The  superscript  '(B)'  denotes  third  order  stress  and  heat  flux  terms  in  the  BGK-Bumett 
equations.  The  coefiBcients  of  third  order  terms  in  the  BGK-Bumett  equations,  fl/s,  are 
given  in  Table  4. 


Y=  1.4 

Y=  1.666 

e, 

2.56 

1.778 

67 

1.36 

1.111 

63 

0.56 

-0.222 

-0.64 

-0.889 

0.96 

0.444 

e« 

1.6 

1.333 

67 

-0.4 

-0.667 

% 

-0.24 

-0.222 

69 

1.024 

1.185 

e,n 

-0.256 

-0.593 

0.) 

1.152 

1.778 

6,7 

0.16 

0.444 

10 


0n 

2.24 

2.222 

0,4 

-0.56 

-1.111 

0,, 

3.6 

3.333 

0,fi 

0.6 

0.333 

0,7 

1.4 

1.667 

0,« 

4.9 

4.167 

0„ 

7.04 

6.222 

-0.16 

-1.778 

07, 

-1.76 

-3.111 

077 

4.24 

4.222 

077 

3.8 

4.333 

074 

3.4 

3.667 

Table  4.  The  coeflBcients  in  the  third  order 
terms  of  BGK-Bumett  Eqns. 


Finally,  the  governing  equations  are  nondimensionalized  by  reference  length  and  free 
stream  variables  and  coordinate  transformed  to  the  computational  ^-t|  domain  by  the 
following  relations 

T  =  /,  * 

^  =  Ux,y), 

and 

Tl  =  Ti(x,;^).  (24) 
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3.  Numerical  Method 

An  explicit  finite  difference  scheme  has  been  employed  to  solve  the  governing 
equations.  The  Steger-Wanning  flux-vector  splitting  meth^<*^  is  appUed  to  the  inviscid 
flux  terms.  The  second-order  central  differencing  scheme  is  applied  to  the  stress  tensor 
and  heat  flux  terms.  In  the  blunt-body  flowfield  calculations  reported  in  this  paper,  fi'ee 
stream  conditions  were  used  along  the  outer  boundary.  First-ordw  extrapolation  of  the 
interior  data  was  used  to  determine  the  flow  properties  along  the  exit  boundary.  Symmetry 
boundary  conditions  were  applied  to  the  stagnation  streamline.  The  first-order  Maxwell- 
Smoluchowski  slip  boundary  conditions®  were  used  on  the  wall  surface  boundary.  The 
first-order  MaxweU-Smoluchowski  slip  boundary  conditions  in  Cartesian  coordinates  are: 


(25) 


(26) 

where 


2-0  / 3  n 


dyj^  ^  A  pT 


.axj  ’ 


and 


2-a  2y  /far] 

'T  —T  _ _ ! - - 

*  ^  a  y  +  lPrv5^j^’ 


The  subscript  ‘s’  denotes  the  flow  variables  on  the  solid  surface  of  the  body.  The 
reflection  coefiBcient,  a,  and  the  accommodation  coefiBcient,  a,  were  assumed  as  1  (for 
complete  accommodation)  in  this  study. 
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4.  Application  to  Blunt  Body 


The  augmented  Burnett  and  the  BGK-Bumett  equations  are  applied  to  compute  the 
hypersonic  flow  over  a  cylindrical  leading  edge  with  nose  radii  of  2m,  0.2m,  and  0.02  m. 
Since  the  numbers  of  grid  lines  are  fixed  in  ^  and  tj  directions,  the  smaller  cylinder  has  the 
finer  grid  system  in  the  physical  domain.  The  grid  system  in  the  physical  domain  is  shown 
in  Fig.  1.  The  flow  conditions  are: 

A/.  =  10.0 
Re.  = 167S 
P.  =23881 
r  =  208.4  ’K 

GO 

r„= 1000.0  'K. 

The  coefficient  of  viscosity  is  calculated  by  the  Sutherland's  law. 

It  =  ^1  - • 

^+^2  (27) 

Various  constants  used  in  the  calculation  for  air  are, 

7  =  1-4, 

Pr  =  0.72, 

P  =  287.04  m*/(scc^ ’A:), 
c,  =  L458  X  10‘  kg  /  sec-  m-’K'^ , 
and 

Cj  =110.4  X 

With  the  given  flow  conditions  and  constants,  the  computations  were  performed  at 
Knudsen  numbers  of  0.001,  0.01,  and  0.1  corresponding  to  the  cylinder  radii  of  2ot,  0.2m, 
and  0.02m  respectively. 
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5.  Results  and  Discussion 


5.1.  Case  1:  (Kn  =  0.1) 

The  comparisons  of  density,  velocity,  and  temperature  changes  along  the  stagnation 
streamline  are  shown  in  Figs.  2,  3,  and  4  respectively.  The  results  are  generally  consistent 
with  those  of  Zhong®.  The  temperature  curves  (Fig.  4)  show  that  the  present  augmented 
Burnett  solution  using  the  Steger-Warming  scheme  has  a  slightly  higher  maximum 
temperature  than  reported  by  Zhong®.  The  Navier-Stokes  solutions  are  also  compared 
with  the  augmented  Burnett  solutions  in  Figs.  2-4.  Since  the  flow  is  in  the  continuum- 
transition  regime  in  this  case,  the  differences  between  the  Navier-Stokes  and  the  Burnett 
solutions  are  significant.  The  shock  width  in  the  augmented  Burnett  solution  is  larger  and 
the  shock  is  upstream  of  that  in  the  Navier-Stokes  solution.  The  density  and  temperature 
contours  of  the  Navier-Stokes  solutions  and  the  augmented  Burnett  solutions  using  the 
present  scheme  are  shown  in  Figs.  5-8.  The  shock  structure  of  the  present  augmented 
Burnett  solutions  agrees  well  vith  that  of  Zhong®.  The  BGK-Bumett  solutions  are 
currently  in  progress  for  this  case. 

5.2.  Case  2:  {Kn  =  0.01) 

The  comparisons  of  density,  velocity,  and  temperature  changes  along  the  stagnation 
streamline  between  the  Navier-Stokes,  the  augmented  Burnett,  and  the  BGK-Bumett 
solutions  are  shown  in  Figs.  9,  10,  and  11  respectively.  The  resulting  curves  are  almost 
coincident  with  each  other.  Only  small  differences  are  observed  at  the  front  of  the  shock. 
The  velocity  curve  of  the  BGK-Bumett  solution  (Fig.  10)  shows  an  unexpected  high  peak 
at  the  front  of  the  shock.  The  density  and  temperature  contours  of  each  equation  solution 
are  also  shown  in  Figs.  12-17.  The  shock  structures  are  also  similar  to  each  other. 

53.  Case  3;  {Kn  =  0.001) 

At  this  small  Knudsen  number,  the  solutions  of  the  Na\aer-Stokes,  the  augmented 
Burnett,  and  the  BGK-Bumett  equations  are  identical.  Since  the  flow  is  in  the  continuum 
regime,  the  Navier-Stokes  equations  already  describe  the  flow  field  accurately.  Figs.  18, 
19,  and  20  show  the  density,  velocity  and  temperature  changes  along  the  stagnation 
streamline  respectively.  Figs.  21-26  show  the  density  and  temperature  contours  for  the 
Navier-Stokes,  the  augmented  Burnett,  and  the  BGK-Bumett  equations. 
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6.  Conclusions 


The  2-D  augmented  Burnett  equations  and  the  BGK-Bumett  equations  have  been 
applied  to  compute  the  hypersonic  blunt  body  flow  (for  at  ^  =  0.1,  0.01,  and  0.001. 
The  explicit  fimte  difference  scheme  with  Steger-Warming  flux-vector  splitting  has  been 
employed  to  discretize  the  convective  terms  in  the  flow  equations.  Simple  second-order 
central  differencing  is  used  to  discretize  the  stress  and  heat-flux  terms.  The  density, 
velocity,  and  temperature  changes  along  the  stagnation  streamline  were  compared  for  each 
set  of  equations.  At  Ah  =  0.1,  the  resulting  flow  properties  and  the  shock  structure  are 
consistent  with  the  results  reported  by  Zhong®.  At  low  Knudsen  number  (Ah  <  0.01),  the 
Navier-Stokes  solutions  and  the  two  Burnett  solutions  are  identical,  "ae  augment^ 
Burnett  equations  were  always  stable  at  all  Knudsen  numbers  and  all  grid  sizes  reported  in 
this  paper.  However,  the  BGK-Bumett  equations  have  experienced  some  convergence 
problem  on  the  finer  grids  at  An=0. 1 .  This  issue  is  being  investigated  currently. 
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Fig.  9.  Density  along  stagnation  streamline  for 
case  2  (^n=0.01). 


Fig.  11.  Temperature  along  stagnation  streamline 
for  case  2  (A!h=0.01). 


Fig.  12.  Navier-Stokes  density  contours  for  case  2 
(X>7=0.01). 
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Fig.  17.  BGK-  Bumett  temperature  contours  for 
case  2  (^=0.01). 


Fig.  18.  Density  along  stagnation  streamline  for 
case  3  (^=0.001). 
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